#include "math.h"
#include <iomanip>
#include <iostream>

double Tik(double (*f)(double), double a, double b, int n)
{
    double x, y;
    double h = (b - a) / n;
    y = 0.5 * (f(a) + f(b));
    for (int i = 1; i < n; i++)
    {
        x = i * h;
        y = y + f(x);
    }
    std::cout << "Tik result = " << std::setprecision(9) << h * y << std::endl;
    return h * y;
}

double Sim(double (*f)(double), double a, double b, int n)
{
    double x, y = 0.0;
    double h = (b - a) / n;
    y = 1.0 / 6.0 * f(a);
    for (int i = 1; i < n + 1; i++)
    {
        x = (i - 1.0 / 2.0) * h;
        y = y + 4.0 / 6.0 * f(x);
        x = i * h;
        y = y + 2.0 / 6.0 * f(x);
    }
    y = y - 1.0 / 6.0 * f(b);
    std::cout << "Sim result = " << std::setprecision(9) << h * y << std::endl;
    return h * y;
}

double Sim(double (*f)(double), double a, double b, double eps)
{
    double y = 0.0, deltay = 0.0;
    int n = 2, maxn = 1000, i = 0;
    y = Sim(f, a, b, n);
    deltay = 10 * y;
    while (deltay > eps && i < maxn)
    {
        deltay = y;
        n = 2 * n;
        y = Sim(f, a, b, n);
        deltay = y - deltay;
        i++;
    }
    std::cout<<"n iter = "<<i<<std::endl;
    return y;
}

double Sim(double (*f)(double), double a, double eps)
{
    double y=0.0,deltay=0.0;
    double b=a,c=10*abs(a+1);
    int maxn=1000, i=0;
    y = Sim(f,b,c,eps);
    b = c;
    c= 10 * abs(b+1);
    deltay = Sim(f,b,c,eps); 
    while (deltay > 0.01*eps && i < maxn)
    {
        y=y+deltay;
        b = c;
        c= 10 * abs(b+1);
        deltay = Sim(f,b,c,eps);
        i++;
    }

    std::cout<<"int a to inf  = "<<y+deltay<<std::endl;
    return y+deltay;
    
}

double intf(double x)
{
    return exp(-x);
}

int main(int argc, char const *argv[])
{
    double y;

    //梯形积分
    Tik(sin, 0.0, 1.0, 100);

    //辛普森积分
    Sim(sin, 0.0, 1.0, 100);

    Sim(intf,0.0,0.0001);

    return 0;
}
